DYNAMIC STRESS ANALYSIS IN THE ROOT REGION 
OF A TURBINE BLADE 


by 

JAY PRASHANT SHRIVASTAVA 


Me 
1*183 
i M 

I 

sh r 

DyN 



DEPARTMENT OF MECHANICAL ENGINEERING 

INDIAN INSTITUTE OF TECHNOLOGY;’ KANPUR 

SEPTEMBER, 1989 



DYNAMIC STRESS ANALYSIS IN THE ROOT REGION 
OF A TURBINE BLADE 


A Thesis Submitted 

in Partial Fulfilment of the Requirements 
for the Degree of 

MASTER OF TECHNOLOGY 


by 

JAY PRASHANT SHRiVASTAVA 


to the 

DEPARTMENT OF MECHANICAL ENGINEERING 

INDIAN INSTITUTE OF TECHNOLOGY, KANPUR 

SEPTEMBER, 1989 



/' 0 1 • M 



CERTIFICATE 


U 


This is 1.0 certify that the thesis entitled, 
’Dynamic Stress Analysis in The Root Region of A Turbine Blade’ by 
Jay Prashant Shrivastava is a bonafide record of work done by him 
under our guidance and supervision and has not been submitted 
elsewhere for the award of a degree. 



(Dr. P. M. DIXIT) 


(Dr. N. S. VYAS) 


Assistant Professor 
Dept, of Mech. Engg. 

Indian Institute of Technology 
Kanpur - 208016 


Assistant Professor 
Dept, of Mech. Ingg. 

Indian Institute of Technology 
Kanpur - 208016 


ii 



ACKNOWLEDGEMENTS 


I express my gra'bi'tude "bo Dr. P. M. Dixit and Dr. N. S. 
Vyas for their inspiring guidance , invaluable suggestions , 
constructive criticism and being constant source of 
encouragement throughout this work. 

I thank Dr. N. N. Kishore for allowing me to use the 
printer of ESA Lab. 

I thank all my friends for helping me throughout my stay 
over here. 

I greatly appreciate the excellent and neat tracings of 
Mr. B, K. Jain. 


— J. P. Shrivastava. 



CONTENTS 


Page 

LIST OF TABLES vi 

LIST OF FIGURES vii 

NOMENCLATURE ix 

ABSTRACT xii 

CHAPTER I INTRODUCTION 1 

1 . 1 In-broduction 1 

1.2 Li^terature Survey 3 

1.2.1 Blade Vibrations 3 

1.2.2 Contact Problem 4 

1.3 Objective and scope of present work 6 

1.4 Plan of the thesis 9 

CHAPTER II THE PROBLEM 10 

2.1 The blade 12 

2.2 The disc 18 

2 . 3 FEM Formulation 20 

2.3.1 Idealization of the contact 20 
problem as a plane stress problem 

2.3.2 FEM Solution procedure 21 

2.3.3 FEM Formulation for plane stress 22 

2.3.4 Formulation of contact problem 27 

2.3.4. 1 Contact conditions 27 

2. 3. 4. 2 Application of contact 31 

conditions 

2. 3. 4. 3 Substructuring 33 

2 . 3 . 4 . 4 Contact conditions 35 

application procedure 

2. 3.4. 5 Evaluation of stresses 36 


iv 



CHAPTER III RESULTS AND DISCUSSIONS 38 

3.1 Stresses in the blade 40 

3.2 Stresses in the disc 44 

3.3 FEM Program for the contact problem 44 

3.3.1 Test problem for the FEM Program 46 

3.3.2 Convergence of the FEM Program 49 

3.3.3 Results for contact displacements 49 
and stresses 

3. 3. 3.1 Mean contact displacements 49 
and stresses 

3 . 3 . 3 . 2 Alternating contact 58 

displacements and stresses 

3. 3. 3. 3 Parametric study 60 

CHAPTER IV CONCLUSIONS AND SUGGESTIONS FOR FUTURE WORK 74 
REFERENCES 76 


V 



LIST OF TABLES 


TABLE TITLE EMK 

3.1 Finite element discretisations for convergence 50 
study 

3.2 Displacement values for different finite 50 


element discretisations 



LIST OF FIGDRES 


FIGURE TITLE PAGE 

NUMBER 


2.1 

2.2 a 
b 

2.3 

2.4 

2.5 

2.6 

3.1 

3.2 

3.3 

3.4 

3.5 

3.6 

3.7 

3.8 

3.9 

3.10 

3.11 


Root domain for the contact problem 11 
Blade mounted on a rotating disc 13 
Bending of a blade element 13 
Rotating disc 13 
Body in plane stress 23 
Discretisation of the domain 23 
Two bodies A and B in contact 28 
The Tee joint between the blade and the disc 39 


Variation of natural frequencies with rotor 41 
speed 

Campbell-diagram 41 

Alternating root stress Vs. rotor speed 42 

Variation of alternating stress along the 43 
blade length 

Variation of mean root stress with rotor speed 43 
Variation of radial stress along the radial 45 
coordinate of the disc 

Variation of circumferential stress along the 45 
radial coordinate of the disc 

Finite element mesh and the deformed shape for 47 
the test problem 

Variation of normal stress along the semi 48 
contact length for the test problem 

Variation of displacements and stresses along 51-52 
the semi contact length for convergence study 


vii 



3.12 Finite element mesh for the contact problem 53 

3.13 Variation of displacements and stresses along 54-55 

the semi contact length for mean stress 

boundary condition 

3.14 Variation of displacements and stresses along 56-57 

the semi contact length for alternating stress 
boundary condition 

3.15 Variation of normal displacement along the 61-62 

semi contact length for different clearances 

3.16 Variation of tangential displacement along the 63-64 
semi contact length for different clearances 

3.17 Variation of normal stress along the semi 65-66 

contact length for different clearances 

3.18 Variation of shear stress along the semi 67-68 

contact length for different clearances 

3.19 Variation of displacements and stresses with 68-70 

X clearance 

3.20 Variation of displacement and stresses with 71-72 

Z clearance 


viii 



NOMENCLATURE 


Area of Cross-section 

Distance from the Neutral Axis 

Proportional Damping Coefficient 
constants 

Young’s Modulus 
Centripetal force 
Shape function 
Forcing function 
Axial thickness 

Moment of inertia 

Blade length 

Number of nozzles 

Natural frequency in i mode 

Time function 

Outer radius of the disc 

Inner radius of the disc 
polar coordinates 
Time 

Total Kinetic energy 

Kinetic energy due to inertia effects 

Kinetic energy due to centripetal effects 

Displacements in X and Z directions 

Strain energy 

Radial displacement 

Potential energy 

Energy dissipated 

Work done by external force 

Cartesian coordinates 

Bending displacement in Y direction 

Blade length considered in root domain 

Distance along the blade length 

Z / 1 

Angle made by tangential direction at the 




[C3 

[D3 

[K3 

[M3 

[N3 

CD3 

{f} 

{F} 

m 

{U} 


contact node with X-axis 

Radial and Circumferential strain 
Density- 

Angular frequency 
Coefficient of friction 
Clearance 

Domain of the problem 

System boundary 

Poisson’s ratio 

Nozzle passing frequency 

Damping ratio 

Angle of pretwist 

Generalised coordinates 

Alternating and mean stress 

Radial and circumferential stress 

Matrix containing derivatives of the shape 

functions 

Damping matrix 

Matrix containing material constants 
Stiffness matrix 
Mass matrix 

Matrix Containing shape functions 
Modal matrix 
Body force vector 

Traction Force Vector 
Force vector 
Forcing vector 
displacement vector 
Strain vector 
Stress vector 


X 



A 

B 



Body A 
Body B 

e Elemental 

i , j , k Dummy indices 

n,s Normal and tangential direction 

T Transpose 


xi 



ABSTRACT 


To determine the dynamic stresses and displacements in 
the root region of a turbine blade, it is not enough to carry out 
the dynamic analysis of the blade alone. Because of the change in 
the geometry and interaction with the disc, this analysis is not 
valid in the root region. Instead, the region consisting of 
appropriate portions of the blade and the disc has to be analysed 
as a contact problem. The individual analysis of the blade and 
the disc provide the stress boundary conditions for the contact 
problem. 

In this work a general computer program is developed to 
analyse first 

(i) the blade as a cantilever beam and the disc as a body in 
axisymmetric plane stress and then 

(ii) a thin layer of the root region as a plane stress problem by 
the finite element method. 

The finite element segment of the program is tested by 
solving a simple problem with known results. This segment is also 
tested for convergence with respect to mesh size. The results are 
presented in the form of normal and tangential displacements and 
normal and shear stresses at the contact interface both for the 
mean as well as alternating stress boundary conditions. The 
variation of the contact displacements and stresses is obtained 
with respect to the following two parameters : 

(a) rotor speed and 

(b) domain geometry. 
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CHAPTER ! 


INTRODUCTION 


1 . 1 Introduction : 

It is widely recognized that fatigue failure of blades is 
a major cause of downtime in turbomachines. Such failures can be 
very costly in terms of lost production and replacement of parts. 
Dewey and Rieger [1] made a survey which revealed that high cycle 
fatigue alone is responsible for atleast 40 % failures in high 
pressure stages of steam turbines. Blade failures due to fatigue 
are predominantly vibration related and suppression of vibration 
must be a high priority item for reducing the number of industrial 
blade failures. 

Blade vibrations are known to be caused by a number of 
different excitation mechanisms such as nozzle passing frequency 
excitation, partial admission operation, stall flutter, etc. Any 
of these sources may cause resonance with some natural frequency 
and give rise to large dynamic stresses. Blade fatigue usually 
initiates in a region of high stress at some structural or 
metallurgical discontinuity. If such conditions sustain for a 
sufficient number of load cycles,- the local fatigue crack may 
continue to grow until failure occurs. 

The earlier design practice has been to avoid fatigue by 
tuning the blade to operate away from the natural frequencies. 
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However in a "typical turbomachine this may not always be possible 
for it may have three to four rotors, each having several stages 
and thus the turbine may altogether have a thousand blades with 
different characteristics. The dynamic stress analysis of the 
blade, then, becomes imperative. The blades in a 
turbomachine are connected to the disc by an attachment region 
called as blade root. The dynamic stress field in blades is known 
to be a function of a number of geometrical and operational 
parameters viz. cross-sectional areas, moments of inertia, 
operational speed, vibrational modes etc. 

In steam turbines the disc being rigid relative to the 
blades, the blade vibration resembles that of a cantilever and 
though the dynamic stresses are dependent on the various geometric 
and operational parameters, it is often likely that the maximum 
stresses would occur at the blade root, especially in the 
fundamental vibratory mode which is most dangerous. Blade root is 
thus, a region vulnerable to failure. The complex shape of the 
root e.g. Tee, Fir-tree, Dove tail causes further stress 
concentration at the bends or sharp corners in the region. 
Moreover, the interfacial friction with the disc in the root 
region and the resulting wear aggravates the vulnerability of the 
root region. 

The blade centrifugal load transmitted to the disc via 
interlocking hooks on either side of the root. The blade root also 
serves the purpose of providing damping, particularly at low 
centrifugal loads, to suppress the critical stresses, both during 
the machine starting/ stopping as well as steady rotor operations. 

All this necessitates stress analysis of the root region. 
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Such a stress analysis would further facilitate comparative 
studies of different types of roots. 


1 . 2 Literature Survey ; 

1.2.1 Blade vibration ; 

Two distinct approaches towards the blade vibration 
problem have been (i) continuum modal approach and (ii) discrete 
modal approach. In the first approach the Lagrange’s equation is 
formulated using the plate or beam theory and the solution of 
equation is obtained directly by approximate methods like 
Rayleigh-Ritz, Collocation, Ritz-Galerkin or by forming a set of 
first order differential equations by linear transformation and 
solving them by numerical integration, as in [23. In the second 
approach, solution is obtained by any one of the following methods 
: Holzer-Myklestad Method, Finite Difference Method, Finite 
Element Method (special case of Rayleigh-Ritz or Galerkin) , 
Transfer Matrix Method , Polynomial Frequency Equation Method etc, 
as in [33- 

While free vibration analysis procedures are fairly well 
established for beam type of blades more work is necessary to 
establish the dynamic behavior of rotating pretwisted plate type 
of blades. A summary paper by Leissa et al [4] compares the 
existing methods for blades with small aspect ratio. 

Blade damping is an important parameter in the fatigue 
design of blades since it is the damping that dictates the 
resonant stress levels. The main damping mechanisms are (i) 
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interfacial damping at the root, (ii) material damping and (iii) 
steam/gas dynamic damping. The superposition of these mechanisms 
determines the stress level achieved in resonant vibration and 
establishes the extent of susceptibility of the blade for fatigue 
failure. Possible contribution of these mechanisms are infinitely 
variable and no predictive theory exists presently to relate these 
damping mechanisms to specific blade damping properties. 

Rao et al [5] devised an experimental method to obtain 
the overall damping coefficients of rotating blades in vibrational 
modes of various rotating speeds and vibrational amplitudes. It 
shows a decrease in damping coefficient with an increase in the 
speed of rotation, the existence of the threshold speed up to 
which frictional damping effect is most predominant and an 
increase in damping coefficients for higher modes of vibration, 
see Rao and Vyas [6] . 

Free vibration characteristics along with damping are 
used to determine the blade response due to a specific excitation 
force. Rieger and Nowak [7] used ANSYS program to determine the 
dynamic stresses. Matsura [8] and Hoyniak and Fleeter [9] have 
done work on blade resonant vibrations. 

1.2.2 Contact Problems : 

Determination of stresses and displacements in the root 
region of the blade is to be treated as a contact problem. For 
simple geometries, analytical methods for contact problems are 
discussed in Gladwell [103 and Muskhelishvili [11]. When the 
geometry of the contact region is complex , it becomes necessary 
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to use the finite element method. Chan and Tuba [121 proposed a 
finite element stiffness method for solving contact problems with 
friction. S. Ohte [13] proposed a method in which the contact 
interface is classified into adhering state and sliding state 
according to the reaction force acting on the contact surface and 
the friction coefficient. Francavilla and Zienkiewicz [14] solved 
the frictionless case of contact problems as a quasi-linear 
problem. N. Okamoto and M. Nakazawa [15] proposed a finite element 
incremental contact analysis with various frictionless conditions, 
considering the contact conditions as additional conditions 
independent of the stiffness equations. Huge and Sauxe [16] used 
FEM and mathematical programming technique to solve frictionless 
contact problems. Here appropriate linearization of the contact 
conditions is adopted and the problem is solved as a non iterative 
optimization problem. Sachdeva and Ramkrishnan extended the method 
of Francavilla and Zienkiewicz to take into account the force 
boundary conditions [17 3 and the frictional effects under 
proportionate loading [18]. Torstenfelt [19] proposed a general 
purpose finite element computer program to solve contact problem 
with friction. Rahman et al [20] used finite element incremental 
procedure for contact analysis where the contact condition at a 
node is established when it numerically penetrates or touches a 
rigid body. Hence one body is considered rigid here. 
Chandrasekhara et al [21] solved the contact problem by imposing 
geometric constraints on the pseudo equilibrium configuration 
defined as a configuration at which the compatibility conditions 
are violated. 
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1 . 3 Ob.iecbive and scope of present, work ; 


The present work is an attempt to determine the dynamic 
stresses and displacements in the root region of a turbine blade. 
The study is restricted to the idealised case of a uniform 
cross-section long blade mounted on a uniform thickness disc 
rotating at a constant angular velocity. By carrying out the 
dynamic analysis of blade as a cantilever beam mounted on a 
rotating support, the stresses can be determined to a reasonable 
accuracy but only in the region sufficiently far away from the 
root. Because of the change in the geometry and interaction with 
the disc, this analysis is not valid in the root region. 
Therefore, the root region is analysed as a contact problem. The 
domain for this problem is chosen to include appropriate portions 
of the blade and the disc surrounding the actual contact area. 
Just like on any other boundary, on the boundary which separates 
the domain from the rest of the blade, either the displacements or 
the stresses have to be known. Here the stresses are used as the 
boundary condition . These stresses are obtained from the dynamic 
analysis of the blade as a cantilever beam mounted on a rotating 
support. For this, it is necessary that this boundary be 
sufficiently far away from the actual contact area, so that the 
Individual blade analysis becomes valid there. The boundary which 
separates the domain from the rest of the disc is also chosen 
similarly and the stresses on this boundary are obtained from the 
analysis of the disc as a body in the state of axisymmetric plane 
stress . 

The domain of the contact problem is a plate-like region 
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of variable -thickBess subjected "to tbe stresses acting in the 
plane of the plate. However, this is not exactly a plane stress 
problem because of the following reason. The dynamic stresses in 
the blade consist of a mean and an alternating part and the 
alternating part varies across the thickness of blade. As stated 
earlier, these stresses become the bovmdary conditions for the 
analysis of the contact problem. Since the stresses on a part of 
the boundary are constant, the problem does not exactly fall into 
the category of the plane stress. However, every small layer of 
the domain can be considered as a plane stress problem. Since the 
alternating stresses are maximum at the extreme layers, one of 
these layers is chosen for the analysis. The chosen layer has 
tensile stresses on the boundary. In summary, a thin layer of 
appropriate portions of the blade and disc surrounding the actual 
contact surface is used as a domain for the analysis of the 
contact problem as a plane stress problem. 

In the blade analysis, the beam theory is used to obtain 
the kinetic and potential energy terms in the Lagrange’s equation 
of motion for free vibrations. The damping is assumed to be 
constant for all modes of vibration. The excitation (nozzle) 
forces are assumed to be harmonic and concentrated at the blade 
tip. The modal analysis is used to obtain the equations of forced 
vibrations which are solved using standard techniques. In the 
disc analysis, the equation of motion is expressed in terms of the 
radial displacements and solved analytically. Because of complex 
geometry, the contact problem is solved by finite element method. 
This method is a modification of the approach suggested by 
Sachdeva and Ramakrishnan [18], They have used the flexibility 
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matrix in their analysis which is obtained by inverting the 
conventional stiffness matrix. However, in the present analysis, 
the stiffness matrix is retained as such. The matrix 

substructuring technique is used to condense the stiffness matrix 
for optimum storage and computational efficiency. Appropriate 
contact conditions involving friction are applied in an iterative 
scheme to the condensed stiffness matrix. Finally the contact 
stresses are obtained by averaging their magnitudes over the 
bodies in contact. 

A general computer program is developed to determine 

(i) natural frequencies, mode shapes, displacements and stress 
fields in the blade, 

(ii) stress field in the disc, and 

(iii) stress and displacement fields in the contact region. 

The blade displacement and stress fields are obtained as 
functions of the rotational speed. The disc stresses are 
represented as functions of the radial coordinates for the various 
rotational speeds. Before applying the finite element segment of 
the program to the contact region, it was tested for a simple 
problem with known results. Further, the convergence of the 
method with respect to the mesh size was checked numerically. In 
the contact domain, the results are presented in the form of 
normal and tangential displacements and normal and shear stresses 
at the contact interface for various values of the rotor speeds. 
The results for mean and alternating boundary conditions are 
presented separately. In addition, the effects of varying the 
contact domain clearances on the interfacial displacements and 
stresses are studied. The variation of clearances is caused 
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either by changing the blade root dimensions or the disc groove 
dimensions . 

This work has several limitations in the form of 
idealizations made and methodologies used. 

1. The blade is idealized as a uniform thickness cantilever beam 
mounted on a rotating support. The actual blade is tapered, 
twisted with an asyimnetric aerofoil cross-section and at times 
may be too short to be treated as a beam. 

2. The disc is considered as a uniform thickness body in the 
state of axisymmetric plane stress. However, because of a finite 
number of blades, it is not exactly axissnranetric. There is also 
some change in the thickness. 

3. The contact region is treated as a plane stress problem. But, 
because of the reasons mentioned above, it is not strictly a plane 
stress problem. It has not been ensured that a small shift in the 
boundaries with the blade or the disc has no appreciable change in 
the contact stresses and displacements. Further, the energy loss 
due to friction is not considered. Finally the iterative scheme 
used for the contact problem, although accounts for friction and 
slip, does not consider the incremental loading. 

1 . 4 Plan of the Thesis ; 

The formulation of the blade, disc and contact problems 
and the solution procedure are described in chapter 2. The results 
and discussion are given in chapter 3, while the conclusions and 
the suggestions for future work in chapter 4. 


8 



CHAPTER H 


THE PROBLEM 


The problem of the dynamic stress analysis of the root 
region of a turbine blade is restricted here to the idealised case 
of a uniform cross-section cantilever beam mounted on a uniform 
thickness disc rotating at a desired constant angular velocity. 
(Fig. 2.1) 

The disc can be analysed as an axisymmetric body 
subjected to a radial centripetal force due to its rotation. The 
rotation of the disc also induces a steady or mean stress on the 
blade. In addition the blade experiences alternating stresses due 
to its vibration on application of lateral external periodic 
forces from the nozzles. The state of stress obtained for such a 
system by analysis of its individual components i.e, the blade and 
the disc would however not hold good in the root region where the 
mutual interaction of the disc and blade and the resulting change 
in geometry have to be accounted for. The root region is hence to 
be treated as a contact problem. The contact domain should 
include appropriate portions of the disc and blade surrounding the 
actual contact area. Such a domain is demarcated by boundaries AB 
and CDEF in Fig. 2.1. The stresses obtained by the individual 
analysis of the blade and the disc are to be used as boundary 
conditions for the contact problem. It is necessary that the 
boundaries of the domain are sufficiently removed from the actual 
contact area for the individual blade and disc analysis to be 
valid at the boundaries. 
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Fig. 2.1 Root domain for the contact problem . 
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2.1 The Blade ; 


The blade is -treated as a uniform thickness cantilever 
beam mounted on a rotating disc Fig. (2. 2a). The kinetic and 
potential energy terms for the blade are setup using beam theory 
and equation of Lagrange "ti employed to generate the eigenvalue 
problem. For ease of analysis damping is assumed to be constant 
for all modes of vibration and the nozzle forces are assumed to be 
harmonic and concentrated at the blade tip. Modal analysis is 
employed to obtain the equations of forced vibrations and the 
stress and displacement fields on the blade are obtained using 
standard techniques. 

The potential energy of the vibrating beam is given by 

1 

V = |eI J (y«)^ dZ (2.1) 

0 

The kinetic energy T, can be seen to be constituted of 
two components - T^, due to inertia effects and T^^, due to the 
centripetal effects . i . e . 

T = Tj + Tq (2.2) 

The component due to the inertia effects readily is 

1 

Tj = I p A J (y)^ dZ (2.3) 

0 

while the component due to the centripetal effects is the work 

done against the centripetal force i.e. 

1 

T^ = J (dF^) u (2.4) 

0 
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where dF^ is -the centripetal force on a small element of thickness 
dZ at a distance Z from the root(Fig 2.2b) 

From Fig 2.2 


Au = dZ (1-cos y’ ) 

dZ „ 

^ ^ (y’ 


Thus 


I 


^ (y’)^ 


dZ 


(2.5) 


0 


The centripetal force on the element is given 


by 


dF^ = p A (R^ + Z)dZ 


( 2 . 6 ) 


Substitution of (2.5) and (2.6) in (2.4) and simplification gives 

1 

T^ = - pio^Aj(y’)^ [R^(l - Z) + -|- (1^- Z^)3 dZ (2.7) 
0 

Assuming a series form displacement field for the blade 


we have 


y = I f.(Z) 


( 2 . 8 ) 


where g.^(t) is a time function and fj^(Z) is the shape function 


given by 


f^(Z) = (i+2)(i+3) 2^"^^ _ i(it3) ^ i(i-H) 


(2.9) 


with Z = 


The assumed shape function satisfies the boundary conditions of a 
cantilever. 

From equations (2.1) to (2.3) and (2.7) to (2.9) we have 

1 


- 4- El I [ 2 ] 


dZ 


( 2 . 10 ) 


0 
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and 


I ^•^|[ 2 Mi] *^2-1 pi^^J^R^d-Z) + 1(1^- 2^)]( 2 


( 2 . 11 ) 


The equation of Lagrange is given by 


A. f ^ 'f 1 

dt ItaJ 


^ T , ^ V _ ^ 

T-^ * -W-^ - ^ 


( 2 . 12 ) 


Substitution of (2.10) and (2.11) in (2.12) yields the eigenvalue 
problem 


[M3{q}+[K3{q}=0 


(2.13) 


where 


M. . = P A f f.t. dZ 
ij 1 j 

0 


; and 


= P A J CRj^d-Z) + -I (1^- Z^)3 f^ f^ dZ 4 ElJ f” f« dZ 


(2.14) 

The eigenvalue problem of equation (2.13) is solved to obtain the 
natural frequencies and mode shapes. 

For ease of formulation damping is ass\imed to 
be viscous and proportional, such that the energy dissipated is 


C T 

W - ^ 

~ p 


(2.15) 


where is proportional damping coefficient. 


The equation of motion for damped vibration can then be 


written as 


[M3{q}+IC3{q>+[K3{q}=0 


(2.16) 
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with [ C 3 = — - — CM] 

The excitation force is periodic with nozzle passing 

frequency v (= ^ 5 ^') • li- is assumed to be harmonic in nature and 

concentrated at the free end of the blade such that 

J'y(2,t) = Fg Sin ut for Z = 1 

= 0 for Z 1 (2.17) 

The virtual work done by this external force is 

1 

<5 = J (Fy <5^) dZ 

0 

= { 00 } Sin ^>t 

■wiiere { } = {F^ (2.18) 

The equation for forced vibration then is 

C M ]{ q‘} + [ C ]{ i } + [ K ]{ q } = { Qg }Sin ut (2.19) 

The forced vibration problem can be solved using modal analysis. 

If { q } = [ 0 ]{ T? } (2.20) 

where [ U ] = Modal matrix 

then 

C M ][ D ]{ >i’} + C C ][ U ]{ ^ } + C K 3C U ]{ ■/? }={Q0}Sin ut 

( 2 - 21 ) 

T 

Premultiplication by [ U ] gives 

C M ]{ 7^'} + C C ]{ 7 ? } + [ K ]{ y? } = { } Sin ut (2.22) 

where 

C M ] = C U 3'^C M ][ D ] 

C C 3 = C U 3'^C C 3C U 3 

C K 3 = [ D 3 ’^C K 3C 0 3 
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(2.23) 


{ } = C U Qg } 

Hence the decoupled equation of motion can be written as 

T? + rj + I? = ^ Sin vt 

for which the solution is 


(2,24) 


where 


'«0N, / «j > 


y> . = 
3 


f 2 , 2.2 . . ^2 „ 2 , 2 

(p -v> ) + 4 ? . p . -u 

^ J J 3 


Sin (r't-V' . ) 


(2.25) 


C. = 




2 «3 ^3 


; and 


. = Tan 


-1 


2 C • p . ‘u 

3 3 


2 „2 

p. - V 


The displacement field can be obtained by using equation 
(2.25) in conjunction with equations (2-20) and (2.8) 

The alternating stress field at any blade cross-section 
at any distance b from the neutral axis of the cross-section can 
be obtained using the formula 

O' = E b y’^ (2.26) 

a 

The mean stress field induced in the blade due to disc 
rotation can be evaluated as fallow . 

The centripetal force on a small element of length dZ at 
a distance Z from the root of the blade is given by (Eq. 2.6) 

dF - p k (R,+ Z) dZ . (2.6) 

C i 


Therefore 


1 

J P A (R^+ Z) dZ 
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= p £o^ A + -| 


2 -,Z=1 


-j 

- 7 .= 


2 = 2 . 


The mean stress experienced by the blade then is 

F 


m A 


2 r 

P R^ (1-Z^) + — 


(2.27) 


2 . 2 The Disc : 


The disc is assumed to be of uniform thickness 
clamped at its centre on a shaft which is rotating at a constant 
angular velocity w about the Y-axis (Fig. 2.3). The disc is in a 
state of axisymimetric as well as plane stress. There is only one 
nontrivial equation of motion. In terms of stress components, it 
is given by 


d 




& 


— + p = 0 


d r 


(2.28) 


Since it is one equation in two unknowns, it 
cannot be solved in the present state. However, it can be solved 
by expressing it in terms of the radial displacements. 

For the axisymmetric plane stress problem, the 
stress-strain relations are 


E 


O' 

r 





E 

^ = g- ( ) (2.29) 

1 - 
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and the strain displacement relations are 


d u 


u 


d r 




(2.30) 


Substituting the relations (2.29) and (2.30) into (2.28) we get 


E 




d r 


d r 


(r 




+ p CO r = 0 


(2.31) 


the solution of this equation is 

2 


u = - 
r 


2 ^2 
p co2 -I- + Cl r -H 


(2.32) 


E 


where Ciand Cg are constants of integration. 

Substituting the expression (2.32) for the radial 

displacement into the equation (2.30) and (2.29) we get the 

following expression for stresses 

E E 3+v ^ ^ 

O' = C. p co*^ (2.33) 

^ 1-v 1+v 8 


E E l+3v ^ » 

^ = C. + — f p co^ r'^ (2.34) 

1-v 1+v 8 

The constants Ci and Cg are to be determined from the 
boundary conditions . Since the shaft to which the disc is 
attached relatively more rigid than the disc, at the boundary i‘=R 2 
the boundary condition becomes 

u^ = 0 at r = Eg (2.35) 

the boundary r = Ri is stress free. Therefore 

0 -^ = 0 at r = Ri (2.36) 

Using these conditions, we get 

1-1^ C„ (3+v)(l-i^) „ „ 

C = 1_ + p R (2.37) 

^ 1+v R 8 E 
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witih 




(l-y^)p 

8E 


(3+v)R^ 


(l+i>)R" 


(l+v)Rf + {1 -i>)R; 


-] 


(2.38) 


Finally, after substituting the values of the constants, the 
expressions for stresses becomes 


a 

r 


(Ei^-r^) [03 


S ^2 

(1+i^) R^^ r^ 


] 




[t 


R^^ + 


E C, 


(1+v) RJ C3 


r2 ^ 


E C. 


(1+1^) C, 


where 


3+y o 

p 

8 


(2-39) 


l+3x> 


3+i> 



(2.40) 


(2.41) 


2 . 3 FEM Formulation : 

2.3.1 Idealization of the Contact Problem as a Plane Stress 
Problem : 

As stated earlier the stresses in the root region of the 
blade are to be determined by treating it as a contact problem. 
The domain of the problem is shown in Fig. 2.1. The domain is a 
plate like region of variable thickness subjected to the stresses 
acting in the plane of the plate. However, this is not exactly a 
plane stress problem because of the following reason. As 

mentioned before, the boundary conditions for the problem come 
from the analysis of the sections 2,1 and 2.2. Specifically, the 
boundary stresses on the boundary AB are obtained from the 
analysis of the blade described in section 2.1, There it is seen 


20 



that, the alternating part of this stress varies across the 
thickness of the domain. Since the stresses on a part of the 
boundary are not constant across the thickness, the problem 
doesn’t exactly fail into the category of plane stress. However, 
every small layer of the domain can be considered as a plane 
stress problem. Since the alternating stresses are maximum at the 
extreme layers, it is appropriate to use either of these layers 
for the analysis. Here the layer with tensile stresses on the 
boundary AB is considered. The analysis of the other extreme 
layer can be carried out in similar fashion but not have done in 
this work. 

Henceforth it should be understood that the domain for 
the contact problem is a thin layer of the blade and the disc in 
the vicinity of the contact surface whose cross-section is shown 
in Fig 2.1. Moreover the stresses on the boundary AB are tensile. 
Because of the complexity of the geometry, the finite element 
method is used for the analysis. 

2.3.2 FEM solution nrocedure : 

In FEM formulation, the governing equations (partial 
differential equations) are converted into an equivalent set of 
integral equations. There are two common procedures available for 
converting the governing equations into equivalent integral 
equations. These are (i) the variational procedure and (ii) the 
Galerkin procedure (or the method of weighted residues). The 
solution domain is divided into several small elements of chosen 
shape and the unknown field variables are expressed in terms of 
assumed approximating functions (or interpolation functions) 



witihin each element. The approximating functions are defined 
within each element in terms of the field variables at specified 
points called as nodes. The nodal values of the field variables 
and the interpolation functions for the element completely define 
the behavior of the field variables within the element. For the 
finite element formulation of the problem, the nodal values of the 
field variables become the new unknown. The integral equations of 
the problem are evaluated within each element using the 
interpolation function for the element. All such elemental 
integrals are assembled to obtain the final matrix equation. The 
unknown nodal variables are calculated by solving the assembled 
matrix equation for the whole solution domain. 

2.3.3 FEM Formulation for plane stress /strain : 

Consider a body O of uniform thickness h bounded by two 
parallel planes Y=-h/2 and Y=h/2 and a closed boundary r (Fig. 
2.4). The part of the boundary on which external stresses act is 
denoted by q . The body is in a state of plane stress . 

In X-Z plane, the strain displacement relations are 

Jin _ ^ . V - — + — 

X '^x ’ 2 “ ■dz ' ^XZ ’’ &X. 

In matrix form 

{ ^ } = [ L 3 { u } (2.42) 

Where { ^ 

{u> = {u w}"^ 


22 




[ L 3 = 


(2.43) 


Then stress strain relationship for a homogeneous, 
isotropic, linearly elastic material under plane stress conditions 


is given by 


Where 


{ ^ } = [ D 3 { ^ } 


{<>■} = {(? T }■ 

'•X z xz 


(2.44) 


[ B 3 = 


1 V 0 

V 1 0 

0 0 ^ 


( 2 . 46 ) 


The strain energy is 


" I ^ 


{£'} {<y}hdA 


(2.46) 






(2.47) 


represent the applied traction on and the body force (per 
unit volume) respectively, then the work done by the external 
forces is given as 


? = f { u } h dA + f { u { f } h dl 


(2.48) 


The expression for Potential Energy of the system is given as 


n = U - W 

P P 


Using Eqs. (2.44), (2.46), and (2.48) 


(2.49) 
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= ^j{ £ }^[ D ]{ £ }hdA - 


I- 


}^{ f }hdl 


(2.50) 


If the domain be divided into N number of elements 

6 


(Fig. 2, 5), over a typical element e let 


{ u } = 


■ 

0 . . 

■ • «n 

0 ■ 

. ® 

N, . . 

. . 0 

N 

n -* 

{ U 

■^nodal 




u. 


w 


1 

e 

e 


u 


w 


"" } 

(2.51) 

be the approximation for the displacement vector, where n is the 

number of nodes per element u^ , > ’ '^n 

displacements at the nodes and are the shape functions 

in the form of polynomials of the coordinates x and z. Then from 
equation (2.42), we get 

{ ^ ) = [ B ] { u } (2.62) 

Where the matrix [ B ] contains the derivatives of the shape 

functions. Substituting Equations (2.51) and (2.52) in to the 

expression (2.50) for the potential energy, we get 
Ne 


n = 

P 


I^J, 

Q 

e=l < 


< - >Ldal t B C D 3 [ B ] t u dA 


N. 


N. 


It 

e=l e 


" >Ldalf « ^ ^ ^ Wal 


[ N 3 { f } 


dl 


(2.53) 
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Here r is the boundary of an element which is common with r 
e q. 

and is the number of such elements. (Fig. 2.5). 

Setting the first variation of n to zero, we get 

Jr 


[ K 3 { U } = { f } 


(2.54) 


Where 


N 


[ K 3 = ^ [ K ]' 
e=l 


(2.55) 


[ K 3®= J C B 3"^ C D 3 [ B 3 h d A, 


a 

f. 

N 


< F > - ^ { F J' 
e=l 


(2.56) 


(2.57) 


{ F }®= J[ N 3’^{ B^ }h dA + I E N 3'^{ f >h dl (2.58) 

D r 

e e 

and { U } is the vector of nodal displacements for the entire 
domain . 

In the equation (2.55) and (2.57), it is understood that 
the sum is performed after expanding C K 3^ and { F to the full 
size. 

The formulation of plane strain problem would be similar 

except that in this case the elasticity matrix C D 3 would be 

1 0 


[ D 3 = 


E(l-i>) 


(1+v) ( l-2i>) 


V 


(1-i^) 

0 0 


1 0 
( 1 - 2 ^^) 


2(l-v) 


(2.59) 
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2.3.4 Formulation of the Gorttact Problem ; 


Two bodies A and B which are in the state of plane 
stress/strain are in contact with each other as shown in Fig 2.6. 
The finite element equations for these bodies are 

['"a] {«a} = {^a} '2.60a) 

[''b] {"b} ={^b} '2-6®=) 

where [ ] . [ Kg ] . { ^ } • ( Cg } , { } , { Fg } are the 

stiffness matrices, displacement vectors and force vectors of 
bodies A and B respectively. The force vector includes the 
contact forces also. 

The equations (2.60) cannot be solved as at the contact 
nodes, neither displacements nor forces are known. To overcome 
this difficulty, contact conditions are used. 

The matrix equation for A and B can be combined into a 
single equation as 


1 "a ] ® 

1 

N 1 


W 1 



1 W j 


w 1 


or [ K ]{ 0 > = { F } (2.61) 

2.3.4. 1 Contact Conditions : 

There are two sets of contact conditions, one on the 
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Fig. 2.6 Two bodies A and B in contact. 
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Let- 


unknown contact forces and the other on the displacements. 

■the node of body A be in contact with the node of body B. 
(Fig. 2.6 ) Let n and s denote respectively the normal and the 
tangential directions at the contact node. Then the conditions on 
forces, as given by Newton’s third law, are 

(^a], = - (^b). '2-62) 

^ •'in ^ jn 

[®'a]. = - (^b], (2.63) 

^ ■'is -*35 

The compatibility of displacement components gives the 

second set of condition 


-'in jn 




(2.64) 

(2.65) 


where £ denotes the clearance 

The compatibility condition for the tangential direction 
is valid only when there is no slip between the contact nodes. 
The slip depends on the magnitude of the normal force at the node 
and the coefficient of friction between the bodies in contact. 
In the case of slip, the condition (Eq. 2.65) on tangential 
displacement should be replaced by the following relation 

1 ( ^b). 1 = H ( ^b). 1 (2 6 

1 V j jsl I V ^ 

This condition can be applied to the body A also. 


29 . 



The contact conditions (2.62)-(2.66) are in ■terms of the 
components with respect to normal and tangential directions, while 
the finite elements equation (2.61) are with respect to 
coordinate system. Therefore, first, these conditions must be 
expressed in terms of the components. From the Fig. 2. 6, we 
get the following relations between the n-s components and 
components . 


= <^A>lx “ * ‘Viz Sin « (2.68) 

Here a is the angle which s-directions makes with x-axis. 
Similarly, the normal and tangential components of other forces 
and displacements can be expressed in terms of the X-z components. 


Substituting the Eq. (2.67)-(2.68) and similar relation for 

"’a'is- 

Eq. (2. 62)-(2. 66) , we get the following conditions 


-‘S'a 

) . Sin a H 

XX 

■ <^A>iz 

Cos a - 

‘^b)jx " 

* ‘Viz 

Cos ct 

= 0 











(2. 

,63) 

(Fa) 

ix 

Cos « + 


Sin a + 

‘Vox 

Cos a 

" <''b)jz 

Sin a 

= 0 











(2. 

.70) 

(Oa) 

ix 

Sin a - 


Cos a - 

‘Db)3x 

Sin a 

■"‘Viz “ = 

£ 











(2. 

.71) 

‘Da) 

ix 

Cos a + 


Sin a - 

‘Db)3x 

Cos a 

- ‘Db)3z 

Sin a 

= 0 



(2.72) 

(ss Cos a + jU sn Sin •+ (ss Sin a - ju sn Cos et)(Fg)j^= 0 

(2.73) 

where ss is sign of (Fg)^^ and sn is sign of (Fg)^^^ 
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Not-e that in case of no slip the condition( 2 . 72) applies while the 
condition (2.73) is to be used when there is slip. 

2 . 3 . 4 . 2 Applic ation of contact conditions : 

It can be seen that corresponding to every contact pair 
(ijj), there are four equations in the set (2.61). Two equations - 
(2i-l)th Eq. for x-direction and (2i)th Eq. for z-direction, for 
the body A and the remaining two (2j-l)th Eq. and ( 23 )th, for the 
body B. At stated earlier, neither the displacement associated 
with these nodes nor the forces (i.e. the right hand side of these 
equations) are known. Therefore, these equations can not be solved 
as such. However for every contact pair, there are four contact 
conditions. Our strategy is to replace these equations with 
unknown right hand sides by the contact conditions which have the 
known right hand side. The matrix operations required for this 
process are described below. 

First we replace the (2i-l)th Eq. of body A by the 
Eq.(2.69). Since (F.).^ and (F. ). are equal to the left hand 
sides of (2i-l)th and (2i)th Eq. of the body A and (B’-n).! and 

• are equal to the left hand side of (2j-l)th and (2j)th Eq. 

Jd 3 

of body B, the condition states that 
(-sinot) (2i-l)th Eq. of body A 
+(cosa) (2i)th Eq. of body A 
+(-5ina) (2j-l)th Eq. of body B 
+(cosa) (2o)th Eq. of body B =0 

Thus after replacement, the (2i-l)th row of the bebomes equal 

to 
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(-sina) ( 2i-l )thi row of [Kj^3 
+ (cosct)(2i) "bh. row of [K^^] 

+(-sina)(2j -1) th row of [Kg] 

+(0030) (2j) th row of [Kg] 

The (2i-l)th row of [Kj^3 > of course, becomes zero. 

Next, the (2i)th Eq. of body A is replaced by the Eq. 
(2.70). The matrix operations for this replacement are similar. 
Now, the (2i) th row of [K^3 is replaced by the linear combination 
of (2i-l)th and (2i)th rows of [K^^] and (2j-l)th and (2o) th rows 
of [Kg], the coefficients of a linear combination given by the Eq. 
(2.70). The (2i)th row of [K^] , of course, becomes zero. Note that 
the Eq. (2.70) involves the original (2i-l)th row of [K^3 and not 
the modified one. Thus (2i-l)th row of [K^] has to be stored 
before applying the first contact condition. 

The third condition is applied in the end. The reasons 
for this change of order will be explained later. So now, 
consider the fourth contact condition. We use this to replace 
(2j-l)th Eq. of body B. In case of no slip the condition is given 
by the Eq. (2.72). Since this Eq. is in terms of displacement 
components, we don’t have to do any row operation like the first 
two conditions while using it for replacement. Instead, we simply 
have to replace the (2j-l)th row of [Kg] by zeros except in the 
columns (2i-l), (2i), (2j-l) and (2j). In these columns, the new 
coefficients will be respectively cosct, sina, (-cosa) and (-sina). 
The (2j-l)th row of [Eg], of course,, becomes zero. IN case of 
slip, the condition is given by the Eq. (2.73). This condition is 
similar to the first two conditions except that it involves only 
two Eq. : (2j-l)th and. (2o) th equations of body B. So here. 
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(2o-l)'th row of [Kg] is replaced by the linear combinations of 
{2j-l)th and (2j) th rows of [Kg], the coefficients of the linear 
combination given by Eq. (2.73). The (2j-l)th row of [Fg] is, of 
course, made zero. Finally, the (2j) th Eq. of body B is replaced 
by the third contact condition (Eq. 2.71). Since it is in terms of 
the displacement components, we simply have to replace the (2j)th 
row of [Kg] by zeros. except in the columns (2i-l), (2i), 

(2j-l)and (2j) where the new coefficients will be sina, (-coset), 

( -sinoi) and coset. The (2j)th row of [Kg] is made equal to & (the 
clearance) . 

If we apply the third Eq. before the fourth, then the 
(2j)th Eq. of body B gets modified. In case of slip, the original 
(2j)th Eq. is needed to apply the fourth condition. So either we 
store the (2j)th row of [Kg] and then apply the third and fourth 
condition in that order or we reverse the order and don’t store 
the (2j)th row of [Kg]. The later alternative is more convenient 
and followed here. 

2. 3. 4. 3 Sub- structuring : 

Since the number of contact nodes is initially unknown, 
an iterative procedure needs to be applied to the Eq. . (2.61). 

Further imposition of contact conditions renders the stiffness 
matrix unsymmetric. An iterative scheme involving large 
unsymmetric matrices requires a large storage and a considerable 
amount of computing time. To reduce both storage and 
computational time, the stiffness matrix is stored in skyline form 
and then decomposed to a much smaller size using substructuring. 
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First step in substructuring is to rearrange the rows and 


columns of the global stiffness matrix as shown below 


>' 1,1 

•^I.II 

^I, III 


^I 



Kn.i 

^II,II 

*^II,III 

I 

Ull 

ZI -1 

*^11 ^ 

*^III,I 

^III,II 

*^III, III 


Uiii 

^ J 


^III 

J 


Here the suffix I denotes the degree of freedom associated with 
contact nodes and the nodes at which non zero forces are 
specified, II denotes the degree of freedom associated with the 
free surface and internal nodes, while III denotes the degree of 
freedom associated with the nodes at which displacements are 

specified. If the prescribed displacements are zero then ^xil“^ 
Further the nodal force at the free surface and internal nodes 
are zero i.e. Fjj = 0 . Thus the Eq. {2.74) reduces to 

%,I ^ = Fj (2.75) 

Kn.I * %,II % = 2 (2.76) 

Eliminating Djj from the above equations we get 

Uj = Fj (2.77) 

where ^II,I ^ (2.78) 

To obtain the matrix Kj j , j is first extracted in 

full form from the global stiffness matrix in skyline form. Next 
Kii is obtained in skyline form from the global stiffness 
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matrix. Then the matrix jj j is obtained columnwise. 

th 

If { C is its i coluum, then it is obtained by solving the 


system 


[ %I,II ] { ^ } ^ { %I.I } 


(2.79) 


i^^Column 


Thus the inverse of jj 

Next the elements of the matrix K, 


is not explicitly found. 


Next the elements of the matrix Kj Kjj ^ are obtained 

by multiplying the rows of jj with the appropriate column of 

I ■ this step as well as in the earlier step, we 

need the matrices Kj jj and j (these are transpose of each 

other). But these matrices are not determined explicitly. 
Instead their rows and columns are generated as and when required. 
Finally Kj j is obtained using the Eq. (2.78). 

Substructuring thus leads to a smaller set of equations 
(2.77), This set is solved after applying the contact conditions. 
Displacements are obtained by back substituting Uj in the 

equation 


^II “ ■ ^II,II ^I 

Once all the displacements are known F^ 

readily found using the third set of equation (2.74). 


(2.80) 


can be 


2 . 3 . 4 . 4 Contact Conditions A'P'plication Procedure ; 


The set of equations (2.77) are solved iteratively by 
applying the contact conditions in the following way 
1. A set of nodal points is assumed to be in contact (it should 
include all possible contact points) and equations (2.77) are 
solved after applying the contact conditions where the fourth 
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condition is the no slip condition (eq. . 2.72). 

2. The normal nodal forces cannot be tensile. Hence all those 
nodes where normal nodal forces come out to be tensile are deleted 
from the possible contact zone in the next iteration. 

3. The ratio of the tangential nodal force to the normal nodal 

force is calculated at all the nodes in contact. If this ratio is 
greater than the coefficient of friction between the two surfaces 
at any node, it implies that this node is slipping in the 

tangential direction relative to the corresponding contact node. 
When a node slips in the next iteration, the contact condition 
(2.72) is replaced by (2.73). The other nodes are still governed 
by the condition (2.72). 

The steps 2 and 3 are repeated until all $.he normal 
forces in the contact zone come out to be compressive and the 

ratio of the tangential force to the normal force for all the 
nodes in contact is either equal or less than the coefficient of 
friction. Once the iterations are over the stresses are computed. 

2. 3. 4. 5 Evaluation of Stresses : 

To evaluate the stress at a point, first the vector 
for the element, to which the point belongs, is obtained 

from the solution vector { U } (the nodal displacement for the 

entire domain). Stresses are calculated by using equation (2.44) 
and (2.52) i.e. 

{<^} = CD3{B}{u (2.81) 

However, we need stresses in the contact region which is 
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not-hing but. a common boundary between the bodies A and B. Our 
approximating function for the displacements is such that the 
stresses do not become continuous across the interelement 
boundaries. Therefore we adopt the following procedure. 

First we choose a set of points in the contact region 
such that the point belongs to only a single element of bodies A 
and B. Then the stress components &, a , r at these points are 

21 !X2lj 

calculated first by assuming that they belongs to the body A and 
then to B. Next the components T and T are calculated by the 

X y 

relation 


T 1 

X 1 

1 

[-1 

X 

T 

X2 


1 (2.82) 

J 

1 

T 

^ XZ 

& 

2 ^ 

1 J 

1 


Here the unit normal vectors for bodies A and B will be 
in the opposite direction. Finally, the normal and tangential 
components of the stress vector { T > are calculated and their 
magnitudes are averaged over the bodies A and B. 


37 



CHAPTER HI 


RESULTS AND DISCUSSION 


A general computer program is developed on the basis of 
the analysis presented in chapter 2. The program takes as inputs 

(i) blade, disc and contact geometry 

(ii) blade and disc material properties 

(iii) rotational speed 

(iv) damping and nozzle excitation data. 

The output consists of 

(i) natural frequencies and mode shapes of the blade 

(ii) stress and displacement fields of the blade 

(iii) stress field of the disc 

(iv) stress and displacement fields in the contact domain. 

The geometry and material properties of the blade disc 
model undertaken for this study is given in Fig. 3,1. 

The blade displacement and stress fields are obtained as 
functions of rotational speed. The disc stresses are represented 
as functions of radial coordinate for various rotational speeds. 
The above stresses at appropriate stations form the boundary 
conditions for the contact domain problem. In the contact domain 
results are presented in the form of normal and tangential 
displacements and normal and shear stresses at contact interface. 
At the interface, points farthest from the central line experience 
maximum normal stresses, which for the present model are most 
critical. As stated in chapter 2 these results are obtained for 
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clearance 


X clearonce 




10 1 I 1^1 5-^ 


Material properties 


Young’s Modulus E =2x10*N/mm^ 
Poisson’s Ratio v =0.25 

Coefficient of friction 
between blade & disc }j=0.2 


Density 

Thickness of the 
blade 


P * 8x10 kg/mm^ 


w = 3 mm 


150 R 


(All dimensions are in mm) 


AB and CDEF boundaries 
define the root domain. 


Fig, 3.1 The Tee joint between the blade and the disc . 



•the ou'bermos'b layer by 'treating the problem as one of plane 
stress. Similar exercise can be carried out for the inner layers 
of the blade-disc assembly. In addition the effect of ■varying the 
contact domain clearances on the interfacial displacements and 
stresses are studied. The variation of clearances is caused 
either by changing the blade root dimensions or the disc groove 
dimensions . 

3.1 Stresses in the Blade ; 

In the program segment for the blade a five term series 

shape function (refer Eq. 2.8) is used which gives converged 

results for free as well as forced vibration response. Fig. 3.2 

depicts the first four natural frequencies of the blade. The 

effect of rotor speed on natural frequencies is also shown. As 

expected the centripetal effect due to rotation causes a 

stiffening of the blade resulting in an increase in natural 

frequency. Fig. 3.3 is the Campbell diagram for the blade. The 

excitation is ass-umed to be a pure harmonic of the Nozzle Passing 

Frequency, NPF (v=n w) and the number of nozzles n is taken as 

s s 

24. The intersection of the NPF line with the natural frequency 
curves gives the resonant rotor speeds which in this case are (602 
RPM, 3285 RPM, 11070 RPM, 23000 RPM) . As can be seen from Fig. 
3.4, which plots the blade root stresses versus rotor speed, the 
intersection of the fundamental vibratory mode with the exciting 
NPF gives the most critical resonance at 602 RPM. The peaks for 
the resonances with higher vibratory modes become successively 
weaker. Fig. 3.5 shows that at 602 RPM the maximvim dynamic stress 
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Fig. 3. 6 Variation of Mean Stress in the Root Region 
with Rotor Speed 





in "the blade occur at t-iie root while at 3285 RPM the maximum 
dynamic stress does not occur at the root but at a point somewhat 
midway between the root and tip. This is so, for the resonance at 
602 RPM is caused due to the interaction between the fundamental 
vibratory mode and the NPP while the resonance at 3285 RPM is 
caused due to the interaction of the second blade mode and the 
NPF. In fact the stress distributions along the blade length for 
various resonant speeds are a consequences of the shape acquired 
by the blade while vibrating in various cantilever modes. 

Fig. 3.6 gives the mean stress experienced by the blade 
due to the rotational effect alone. The centripetal effect being 
more pronounced at higher speeds the mean stresses on the blade 
increase with an increase in rotational speed. 

3. 2 Stresses in the Disc : 

The disc is in a state of axisymmetric plane stress. The 
radial and circumferential stresses as functions of radial 
coordinate are shown in Fig. 3.7 and 3.8 for a typical rotor speed 
of 602 RPM. 

The blade and disc stresses thus obtained for a desired 
rotor speed at the boundaries AB and CDEF shown in Fig. 3.1 are 
input boundary conditions to the program segment for the root 
domain analysis. 

3. 3 FEM Program for the Contact Problem : 

The program segment for the root domain is based on the 
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Radiol Coordinate, mm. 


Fig. 3. 7 Variation of Radio! Stress along the 
Radial Coordinote of the Disc. 



Fig 3-8 Voriotion of circumferential stress along the 
Radial Coordinate of the Disc. 


finite element formulation described in chapter 2 to solve plane 
stress (or plane strain) contact problem with friction. The 
program uses four noded isoparametric elements and gives the 
displacement and forces at each node as well as normal and shear 
stress variation along the contact length. 

3.3.1 Test Problem for the FRM Program : 

Before application of the program to the blade root 
domain it was tested for a simple plane strain problem with known 
results. Two prismatic bodies of rectangular cross-section in 
contact and with a point load at their midspan were considered 
(Fig. 3.9). The finite element mesh for body A consists of 12 
elements and 21 nodes while that of body B consists of 30 elements 
and 44 nodes. To begin with 7 nodes are assumed to be in contact. 
Four iterations are needed to converge to the correct solution. 
The converged solution shows that only 5 nodes are in contact. 
The displacements and stresses at the contact surface are in good 
agreement with those of Ramesh [25] . The plane stress version of 
the problem is solved by Sachdeva and Ramkrishnan [18] and the 
results for the plane stress version with present program match 
quite well with their results. The deformed shape of the two 
bodies in contact (on an enlarged scale) and the normal stress 
variation along the semi contact length are shown in Figs. 3.9 
and 3.10. 
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Load-P .ISSON 



Material properties Deformation Scale 

Young's Modulus E = 0.931 x 10 N/mm^ 0 40 80 120 160 pm 

Poisson’s Ratio = 0.25 
Coefficint of friction p = 0.20 

Fig. 3.9 Finite element mesh and deformed shape for the 
test problem. 
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Fig. 3.10 Variation of normal stress along the semi - 
contact length for the test problem. 
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3.3.2 Ga nvergence of -bhe FEM ; 


Since the contact problem was a non-linear problem it is 
necessary to check the convergence of FEM Program results 
numerically with respect to mesh size. The convergence study was 
made for a rotor speed of 2000 RPM. The program was executed for 
the discretisations listed in Table 3.1 

The variation of the normal and tangential displacements 
and stresses for all three discretisations are shown in Figs. 
3.11. From the curves a very good and fast trend towards 
convergence can be seen. The numerical values of these quantities 
at selected points along the contact length are shown in Table 
3.2. 

The third finite element mesh shown in Fig 3.12 is chosen 
for carrying out further root domain analysis. 

3.3.3 Results for the Contact Dis-placements and Stresses : 

3 . 3 . 3 . 1 Mean Contact Displacements and Stresses : 

The results for stresses and displacements for the 
contact interface of the blade disc assembly are given in Figs. 
3.13 and 3.14. The results are presented for the semi contact 
length marked ba and in the Fig. 3.1. The lengths are measured 
from the Z-axis shown in Fig 3.1. Figures 3.13 depict the normal 
displacement, tangential displacement, normal stress and shear 
stress variation along the semi contact length for various speeds 
of rotation, when the mean stresses from the blade are given as 
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Table .3.1 

Finite element discretisation f,or convergence study 


Discretisation 

Blade part 

Disc part 

Number of 

elements 

Number 

of nodes 

Number of 

elements 

Number 

of nodes 

I 

26 

42 

40 

63 

II 

34 

54 

48 

75 

III 

54 

1 

1 

i 

83 

70 

106 


.Tab le : 3,2 

Displacement values for different finite element discretisation 


Distance 
from the 
centre line 
i . e . z-axis , 
mm. 

Normal 

displacement 

xl0 ^mm. 

Tangential 

7 . xl0 mm. 

displacement 

jSt 

Mesh 

Iind 

Mesh 

Mesh 

jSt 

Mesh 

Iind 

Mesh 

IlfSt 

Mesh 

7.75 

20.89 

21.043 

21.106 

2.0521 

2,0476 

2.0471 

9.00 

17.161 

17.209 

17.225 

3.8852 

3.8880 

3.8933 

10.0 

13.303 

13.256 

13.239 

4.9543 

4.9779 

4.9901 
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Fig.3.11(a) Variation of Normal Displacement along the 

Semi Contact Length for three Discretisations. 
(Convergence Study) 



Fig-3.1 1(b) Variation of Tangential Displacement along 

the Semi Contact Length for three Discretisations. 
(Convergence Study) 
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d-4— Elements in contact 
o- 8-“Elem«nt8 in contact 
«-1 6-Elements in contac 


Semi-contact length, mm. 


Fig.3. 11(c) Variation of Normal Stress along the Semi 
Contact Length for three Discretisations. 
(Convergence Study) 


9 4-Eiements in contact 
o 8-Elements in contact 
« 16-Element5 in contact 


Semi— contact length, mm. 


Fig.3. 11(d) Variation of Shear Stress olong the Semi 
Contact Length for three Discretisations 
(Convergence Study) 







Fig.3.12 Finite element mesh for the contact problem. 
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Semi contact length, mm. 



3 13(a) Variation of Normal Displacement along the 
Contact Length for various Rotor Speeds 
(Mean Stress Boundary Condition) 


Semi 



Scmi contoct length, rnm. 


Fin 3 13fb) Variation of Tangential Displacement along the 
Fig 3. Contact Length for various Rotor Speeds, 

(Mean Stress Boundary Condition) 
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Fig 3J3(c) Variation of Normal Stress along the Sem 
Contact Length for variour Rotor Speeds 
(Mean Stress Boundary Condition) 



Fig 3.13(d) Variation of Shear Stress along the Semi 
Contoct Length for various Rotor Speeds 
(Mean Stress Boundary Condition) 
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11070 R.P.M. 


23000 R.P.M. 



10 T -r n ' Trrn ( \ r - n-n m v i t \ i Ti n"rT’ n ' rn -'r n T f ’f n n r nTT '" | 

7*50 8.00 8.50 9.00 9*50 10.00 

Semi contoct length, mm* 


Fig 3.14(a) Variation of Normal Displacement along the Semi 
Contact Length for various Rotor Speeds, 
(Alternating Stress Boundary Condition) 



Fig 3.14(b) Variation of Tangential Displacement along the Semi 
Contact Length for Various Rotor speeds 
(Alternating Stress Boundary Condition) 
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Normal Stress, N/mmf 


1C 


II 


1 


Semi contact length, mm. 

Fig 3.1 4(c) Variation of Normal Stress along the Semi 
Contact Length for various Rotor Speeds. 
(Alternating Stress Boundary Condition) 



Fig 3.14(d) 


Variation of Shear Stress along the Semi 
Contact Length for various Rotor Speeds 
(Alternating Stress Boundary Condition) 




boundary conditions to the contact domain. It can be seen from 
these figures that the displacements and stresses increase all 
along the semi contact length with an increase in rotational 
speed. This is so because with the mean stress from the blade as 
a boundary condition, only centripetal forces due to rotation are 
being accounted for. Higher centripetal forces result in greater 
displacements and stresses at the contact interface. It can also 
be observed from these figures that at any particular rotor speed, 
the normal displacement of the semi contact interface is higher in 
the regions which are closer to Z-axis as compared to those in 
regions located farther from the Z-axis. On the other hand the 
normal stresses are lower in regions closer to Z-axis as compared 
to those in region farther from the Z-axis. This can be explained 
on the basis of the fact that the interlocking ’hook’ abed, of the 
disc in Fig. 3.1 behaves somewhat like a cantilever experiencing a 
bending load from the blade root, allowing the blade root a 
greater displacement towards its ’free end’ be as compared to its 
fixed end ad and resulting consequently in lower stresses in 
regions towards be as compared to those in regions closer to ad. 

3 • 3 • 3 . 2 Alternating . Contact Displacements and Stresses : 

Similar trends at any particular rotor speed can be 
observed in Figs. 3.14 where the alternating stresses from the 
blade now form a boundary condition for the contact region. 
However in this case, unlike the results presented in Figs. 3.13, 
the stress and displacements values along the semi contact length 
need not necessarily increase with increasing rotor speed. The 
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st.resses at. any point on the semi contact length can be seen to be 
higher for a rotor speed of 602 RPM as compared to those at 3285 
RPM. The RPM of 602 is a resonant rotor speed caused due to the 
interaction of the blade’s first vibratory mode with the nozzle 
passing frequency. The resonance at 3285 RPM is caused due to the 
interaction of second vibratory mode with nozzle passing 
frequency. The alternating bending stresses on the section AB of 
the blade which is at a distance of 20 mm from the root (Fig. 3.1) 
are higher at 602 RPM as compared to those at 3285 RPM (refer Fig. 
3.5). Hence these stresses which form the boundary conditions for 
the root domain bounded by AB and CDEF result in higher stresses 
at contact interface fea for 602 RPM and lower stresses for 3285 


RPM. 


From the beam theory the maximum stress at the root is 
found to be 940 N/mm^ (Fig. 3. 5) while its value is 3025 N/mm when 
the root region is analysed as a contact problem (Fig. 3.14). The 
above difference in stresses clearly emphasizes the need for the 
analysis of the contact region. On the basis of the analysis and 
results presented, it is possible now to obtain the mean and 
alternating stresses at the blade-disc root interface. The 
vulnerable regions of the interface can thus be identified. The 
mean and alternating stresses would form the two axis of the mean 
stress diagram (i.e. Goodman diagram) to predict the fatigue 
damage and residual life of the vulnerable interfacial zones. 
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3.3. 3. 3 Earametrlc St.ndv : 


A study was further carried out to examine the effect of 
variation of X and 2 clearances, at the blade root junction (See 
Fig. 3.1), on the interfacial stresses and displacements. The 
variation in the clearances is caused either by 

(i) keeping the disc dimensions fixed and varying the blade 
dimensions. 

or by 

(ii) keeping the blade dimensions fixed and varying the disc 
dimensions . 

These results in terms of displacements and stresses for 
the entire semi contact length are presented in Fig. sets 3.15 
3.18 for various X and Z clearances. To display the trends more 
clearly the variation of the stresses and displacements at the 
mid-node alone of the semi contact length with variation in X and 
Z clearances are presented in Fig. sets 3.19 and 3.20, The effect 
of once keeping the blade dimension fixed and then keeping the 
disc dimensions fixes is presented simultaneously on the same 
graph. 

Starting with an initial set of X and Z clearances (0.25 
mm. and 0.25 mm. respectively). Figs, 3.19 and 3.20 depict the 
effect of reduction in clearances. It can be seen that a 

reduction in X clearance ’ in general ’ leads to a reduction in the 
contact displacements and normal stresses. This is so for a 
reduction in X clearance implies eibher an increase in blade root 
cross-section area /£> (when disc dimensions are fixed) or an 
increase in contact area t>a shown in Fig. 3.1 (when the blade 
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Normal Displacement, 



Fig 3.15 (q) Variation of Normal Displacement along the Semi Contact 
Length for Various Values of X Cirarences. 

• (Blade Dimensions Fixed). Rotor Speed = 602 RPM. 



Fig. 3. 15(b) Variation of Normal Displacenemt along the Semi 
Contact Length for various values of X Clearances. 
(Disc Dimentions Rxed). Rotor Speed = 602 RPM. 


6l 



Normal displacement, mm. 



Fig.3. 15(c) Variation of Normal Displacement along 
Contact Length for various values of Z 
(Blade Dimensions Fixed). Rotor Speed 



Fiq.3. 15(d) Variation of Normal Displacement along 
' Contact Length for various Values of 

(Disc Dimensions Fixed) Rotor Speed 


the Semi 
Clearances. 
= 602 RPM. 


0.25 mm. 
0.15 mm. 
0.10 mm. 
0.01 mm. 


I the Semi 
Z Clearances. 
= 602 RPM. 
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Tangential Displacement, mm. 


Semi Contact Length, mm. 


Flq.3. 16(a) Variation of Tangential Displacement along the Semi 
Contact Length for various Values of X Clearances. 
(Blade Dimensions Fixed). Rotor Speed = 602 RPM. 



Fiq. 3.1 6(b) Variation of Tangential Displacement along the Semt 
Contact Length for various Values of X Clearances 
(Disc Dimensions Fixed). Rotor Speed = 602 RPM. 





Fiq.3- 16(c) Variation of Tangential Displacement alone 
Contact Length for various Values of ZJ. 
(Blade Dimensions Fixed). Rotor Speed — 


the Semi 
learances. 
602 RPM. 



Fig.3.1 6<S^ariation of Tangential Displacemant along the Semi 
Contact^ Length for various Values of Z Clearances. 
(Disc Dimensions Fixed). Rotor Speed = 602 RPM. 
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Fig 3.17(o) Variation of Norma! Stress along the Semi Contact 

Length for various Values of X Clearances 
(Blade Dimensions Fixed). Rotor Speed = 602 RPM. 



Fin 3 17(b) Variation of Normal Stress aionq ^e Semi Contact 
^ ' Length for various Values of X Clearances 

(Disc Dimensions Fixed). Rotor Speed = 602 RPM. 
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Fia 3 17(c) Variation of Norma! Stress along the Semi Contact 
Length for various Values of Z Clearances 
(Blade Dimensions Fixed). Rotor speed — 602 RPM. 



Fig 3.17(d) Variation of Normal Stress olong the Semi Contact 
Length for various Values of Z Clearances. 

(Disc Dimensions Fixed). Rotor Speed = 602 RPM. 
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Fig.3.1 8(a) Variation of Shear Stress along the Semi Contact 
Length for various Values of X Clearances. 

(Blade Dimensions Rxed). Rotor Speed = 602 RPM. 



(Disc Dimensions Fixed). Rotor bpeea 
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1100 



Fig.3. 18(c) Variation of Shear Stress alonq the Semi Contact 
Length for various Values of l Clearances. 

(Blade Dimensions Fixed). Rotor Rxed = 602 RPM. 



Fig-3.1 8(d) Variation of Shear Stress along the Semi Contact 
Length for various Values of Z Clearances. 

(Disc Dimensions Rxed). Rotor Speed == 602 RPM. 




Fig. 3. 19(a) Variation of Norma! Displacement at the Mid— Node 
of the Semi Contact Length against X Clearance. 



X Clearance, mm. 


Fig.3. 19(b) Variation of Tangential Displacement at the Mid- 

Node of the Semi Contact Length against X Clearance. 
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[•li] 


X Clearance, nnm. 


Fig. 3. 19(c) Variation of Norma! Stress at the Mid— Node 

of the Semi Contact Length against X Clearance 


X Clearance, mm. 


Fig-3.1 9(ci) Variation of Shear Stress at the Mid— Node 

of the Semi Contact Length against X Clearance 






Fig.3.20(a) Variation of Normal Displacement at the Mid-Node 
of the Semi Contact Length Against 2 Clearance. 


E 0.028 - 
E 



Disc dim ensions fixed 
Blade dimensions fixed 


0-024 -j — r— I — 1 1 " I 'I I 1 — r "; "T ' i — i — i " " i i r " i fj- . 

0.00 0.10 0.20 
2 Clearance, mm. 


T" 1 { 

0.30 


Fig. 3. 20(b) Variation of Jangential Displacement at the Mid— 
of the Semi Contact Length against 2 Clearance. 
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Normal Stress 



Z Clearance, mm. 


Fig. 3. 20(c) Venation of Normal Stress at the Mid-Node 

of the Semi Contact Length agoinst 2 Clearance 



Z Clearance, mm. 


Fig.3,20(d) Venation of Shear Stress at the Mid— Node 

of the Semi Contact Length against 2 Clearance. 
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dimensions are fixed). Both result in an increase in load bearing 
capacity consequently causing reduction in the interfacial 
displacements and normal stresses. It can also be seen that the 
increase in load bearing capacity is more prominent when the 
contact area ba. is increased as compared to an increase in the 
blade cross-sectional area at /£». It can be seen however that 
after a certain point an increase in the contact area ha (by 
varying the disc dimensions) results in increased displacement. 
This can be explained on the basis of the fact that an increase in 
contact area ha also implies a reduction of stiffness of the 
cantilevered 'hook' abed.. Stiffness reduction becomes more 
prominent after a certain point whilst increasing the contact area 
for increasing load bearing capacity. 

With the decrease in Z clearance by changing the disc 
groove dimension, the stiffness of the ’.hook’ <2bcdincreases which 
causes decrease in the normal displacement. But the stresses 
remain unchanged. If Z clearance is varied by varying blade root 
dimensions, there is no appreciable change in displacement or 
stresses . 

It is realized that a more detailed parametric analysis 
and consideration of different types of roots (e.g. Firtree, 
Dovetail etc.) is essential to gain more meaningful insights. 

However the present study does underline the relevance of 
such future work, if undertaken. 
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CHAPTER IV 


CONCLUSIONS AND SUGGESTIONS FOR FUTURE WORK 


The general computer program developed in chapter three 
determines 

(i) the natural frequencies, mode shapes and the displacements 
and stresses as functions of the blade coordinate for various 
rotor speeds , 

(ii) the radial and circumferential stresses as functions of the 
radial coordinate for various rotor speeds, and 

(iii) the normal and tangential displacements and the normal and 
shear stresses along the contact surface for various values of the 
following two parameters 

(a) rotor speed, and (b) domain clearances. 

Based on these results following conclusions can be 

drawn . 

(i) From the beam theory, the maximum stresses at the root is 

found to be 940 N/mm while its value is 3025 N/mm. when the root 
region is analysed as a contact problem. There is an increase in 
the maximum root stress by about 300 % . This difference is 

clearly emphasises the need for the analysis of the root region as 
a contact problem. 

(ii) For mean stress boundary condition, the contact displacements 
and stresses increase with the rotor speed. On the other hand, 
for the alternating stress boundary condition, there is no such 
simple trend. Here the variation with respect to rotor speed 
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depends on the dominant mode of vibration. 

(iii) ’In general’, the reduction in X clearance leads to a 
reduction in the contact displacements and the normal stresses. 
This reduction is more prominent when the clearance is reduced by 
keeping the blade dimensions fixed. Comparatively, the change in 
Z clearance has a very little effect on the contact displacements 
and stresses. 

As mentioned in chapter one, this work has several 
limitations in the form of idealizations made and methodologies 
used. In future, some or all of these can be accounted for. In 
particular the following modifications can be undertaken : 

(i) In the analysis of the blade, the actual geometry of the blade 
should be considered. It means the effects of taper, twist, the 
asymmetry of the cross-section and the aspect ratio should be 
accounted for. 


(ii) Because of the finite number of 

blades 

and 

the change 

in 

geometry resulting from the grooves 

, the 

disc 

should 

not 

be 

considered as an axisymmetric body. 

Instead, 

its 

typical 

sector 


should be analysed as a plane stress problem. 

Instead of treating a small layer of the contact region as a 
plane stress problem, the whole of it should be analysed as a 3-D 
Problem. Further the boundary of the domain should be chosen such 
that it has no effects on the contact displacements and stresses. 
The energy loss due to friction also has to be accounted for. 
Finally, for transient problems, the incremental loading procedure 
should be adapted. 
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